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system will survive or complete an intended mission of 
certain duration. 

Describing all possible ways that a system can survive a 
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reliability computations. Reliability computation for a 
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I. INTRODUCTION 



In making' a mathematical model for a real-life phenom- 
enon it is always necessary to make certain simplifying 
assumptions so as to render the mathematics tractable. One 
of the simplifying assumptions that is often made is to 
assume that certain random variables are exponentially 
distributed. The reason for this is that the exponential 
distribution is relatively easy to work with. 

The property of the exponential life distribution which 
makes it easy to analyze is its (memoryless) lack of deteri- 
oration with time. By this we mean that if the life time of 
an item is exponentially distributed, then an item which has 
been in use for a certain amount of time is as good as a new 
item in regards to the amount of time remaining until it 
fails . When the life time of an item is exponentially 
distributed, the failure rate function for the item is 
constant . 

Under the assumption of constant component failure 
rates, it is possible to build a reliability shorthand 
[Ref. 1] for any system. The term system is used to describe 
a set of components organized to perform some mission. 
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Any study on system reliability requires a description 
for the system's life and a derivation of the system's 
survival function. The reliability shorthand gives a simple 
and easy way for describing a system's life, but it is 
difficult to implement computationally since it involves 
considerable complexity in handling convolutions. Here, the 
term convolution refers to the summation of independent 
random variables (lives). 

This paper presents an analytical approach for obtaining 
a general equation for the survival function of any convolu- 
tion of independent and exponentially distributed random 
variables . 

Section 2 deals with convolutions in detail and gives a 
mathematical derivation for a general equation for the 
survival function for any convolution of exponential random 
variables, using Laplace transforms and Theorem of Residues. 

Section 3 gives a mathematical derivation for an alter- 
native formula for the survival function of any convolution 
of exponential random variables by computing the coeffi- 
cients of all the polynomials that accompany the exponential 
terms in the survival function equation. 

Appendix A contains a computer program written in 
Fortran for computing the reliability of any convolution of 
inependent and exponentially distributed random variables. 
This program uses the general equation in Section 2. 
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Appendix B contains another computer program written in 
Fortran to compute the reliability of any convolution of 
exponential random variables, using the general formula in 
Section 3 . 
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II. CONVOLUTIONS 



Under the assumption of constant component failure 
rates, this section presents an analytical approach for 
finding a general equation for the survival function of any 
convolution of independent and exponentially distributed 
random variables. 

i 

A. SURVIVAL FUNCTION BY INTEGRATION 

Let us start with a simple system; a standby system 
having one active component and one cold spare component . 
The life time of the system is 

T = Tj + T£ , 

where 

T, ~ EXP ( X | ) , 

T z ~ EXP ( X 2. ) . 

and 

Tj, T^ are independent. 

An active component A is to complete a mission of dura- 
tion t, while a spare component S replaces the active compo- 
nent when it fails. The life duration of the active 
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TIME 



WORKING COMPONENTS 



Figure 2.1 A Single Active Component With One Spare 

component is T j and the life duration for the spare compo- 
nent, if it is used, is T 2 (see Fig 2.1). 

The system survival function is given by 

F T (t) = P ( T > t) 

t 

= ? T (t) * / F X (t-s).f T (s) ds 



-X, t r -X,s 

e + / e \ e ds , 



where 



f T t ( s ) 



- dF T ^(s)/ds 

X, r X|S . 



By integration, the survival function for the system is 



— / s X 2. — Xjt X, — X^"^ 

F T (t ) = e + — 1 — e , t>0. 



X r X, 



X, X2_ 
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Now, suppose a system consists of one active component 
with constant failure rate \ ^ and two spares having constant 
failure rates X £. anc * X3 res P ect ively . The life time of the 
system is 



where 



and 



T = T, + T2 + T 3 , 

T, ~ EXP(\, ) , 

T z ~ EXP(X Z ) . 

) 

T 3 ~ EXP ( X 3 ) . 

T | , T ^ > T 3 are independent . 



The shorthand notation for the system is 

J I 

EXP ( X | ) + EXP ( X £ ) + EXP(X 3 ) • 

The survival function of the sustem is given by 



Fy (t) - F Tj+-T£( t ) + J F T^( t_s ) • f TjrT^( s ) ds 



where 



Vt> (t) 



F-^(t-s) 




X| -X,t 

• e 4 , 

X r Xz 
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and 



f lj + T z ( s ) = - dF Tj + T z ( s )/ ds 




By integration, the survival function for the system is 
F T (t) = F T| +T z f-^(t) 

XzX 3 e X|t | X,X 3 e-^ t | X.X^e^ t?Q 
(\fXi)(Xr\) (XrXz)(Xj~Xz) (Xf~ Xj)(Xz 'X^) 

To compute the survival function for a system of n 
lives , 

T = T, + . . .+ T n , 

we can proceed in a similar manner. However the integrations 
will become increasingly complex, particularly if the tacit 
assumption made so far that all component failure rates are 
different is abandoned. Therefore, we will use another math- 
ematical way to compute the survival function without direct 
integration. This will be done with Laplace transforms and 
the Theorem of Residues in the next subsections. 
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B. CLASSIFICATION OF CONVOLUTIONS 



Consider a system which has a life time 
T = T, + ... + T n , 

where T, , ..^T^ are independent and exponentially distrib- 
uted random variables. This system can be classified as one 
of three possible convolution cases : 

1. First Case 



The system has n dissimilar 


failure 


rates , 


Xj 4 Xj , i t j , i , j 


= 1, •• 


. , n. 


Second Case 






The system has n identical 


failure 


rates , 


? — 1 
II 

•n 

— > 
II 


. , n. 




Third Case 







The system has a combination of the first case and 
the second case, some similar and some dissimilar failure 
rates . 

We consider the third case to be of particular 
interest, since it better represents the general situation 
in computing the reliability of any convolution of exponen- 
tially distributed random variables. 
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C. CHARACTERISTIC FUNCTIONS AND LAPLACE TRANSFORMS 



Consider a random variable T which has a distribution 
function Fy(t) and a density function f j (t). The charac- 
teristic function of T is defined, for any real number u, 
[Ref. 2] by 




= E 



[e'^] 



+ oo 



/ 



- QO 



iut 

e dF T ( t ) 



+ 00 



iut 



f j (t) d t 



If T is non negative (T >0) and exponentially distributed 
with a failure rate \ , then 

OO 

r f iut 

0 T (u ) - J e • f T (t) dt 

0 

iut v -\t 

e \ e dt 

0 

X 

(X-' u ) 
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Now, if 



T = T j + ... + T n , 

where 

T • — EXP ( X j ) , i = 1 , • • • > n , 

and Tj,...,T n are independent, then the characteristic 
function of the sum (convolution) T^+ ... + is 




This implies that the characteristic function for the sum 
Tj+ ... + is just the product of the individual charac- 

teristic functions. 
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When dealing with random variables which only assume non 
negative values, it is more convenient to use Laplace trans- 
forms rather than characteristic functions. The Laplace 
transform of the random variable having distribution Fy(t) 
is defined [Ref. 2] by 



<^) T (s) = E 



-St 



oO 



- s t 

e dF T ( t ) 



oO 



•/ 



-sf 

e f T (t) dt , 



where f -p (t) is the density function of T. This integral 
exists for a complex variable s=a+bi, where a^O . 

If T, , . . . ,T n are independent random variables, then the 
Laplace transform of the sum (convolution) T ( +...+ T n is 
given by 

r-sC T i + '" +T n) 



Cf> (s) = E 
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n 



] r cb -T- ( s ) . 

i=i *» 



It is important to note that the Laplace transform uniquely 
determines the distribution. That is, in fact, there is a 
one to one correspondence between distribution functions and 
Laplace transforms. 

Now, we will consider the following question. Given a 
function <p . j-(s), does there exist ,a function f-p(t) whose 
Laplace transform is <^)-p(s), and if it exists, how can it 
be determined ? The answer to the first question is not 
always positive. In - general , the function <j6-p( s ) must 
satisfy some restrictions which will be discussed in the 
next subsection. If the answer to ' the first question is 
positive, then the inversion of the Laplace transform by a 
complex integral gives the answer to the second question; 
That is, the density function for the random variable T is 
given by [Ref. 3] 



f T (t) = 




T 



(s) 



c*+ loo 

I r st 



2. TT i 

o( — i oo 



f St , 

I e . -j- ( s ) ds , 



19 



where the notation ( s )] refers to the inversion of the 
Laplace transform and the path of the integral above is the 
straight line Re s=°< (o< >0) parallel to the imaginary axis. 



D. LAPLACE INVERSION AND THE THEOREM OF RESIDUES 



As was stated in last subsection, if we have a Laplace 
transform (^6-j-(s) for a non-negativ$ random variable T, then 
the Laplace inversion T[4> y(s)j gives the density function 
f T (t). 

When dealing with a convolution of random variables, 
that is, T - T, + ... + T n , where T ( , . . . , T n are independent, 

the convolution density function is given by [Ref. 3] 

-I 



f T 00 = X[0T(s>] 



2TT| 



r / 



oi + I <*> 

st A 

e . ©j(s) ds, 



<*— i 



where (p y(s) is the Laplace transform of the sum 
(convolution) T, + ... + T n . 

To obtain the convolution density function, it is clear 
from the above integral that we have to evaluate a complex 
integral, but the Theorem of Residues reduces the evaluation 
of a complex integral along closed contours to passages to 
the limit and differentiations. 
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Suppose a function G(s) is holomorphic in a finite 
domain C, the term holomorphic meaning that the function is 
analytic and single- valued . Then the integral 

— m G(s ) ds 

ZTT i J 

evaluated along it's boundary, vanishes. If the function 
G(s) is not holomorphic in C, the integral may not be zero. 
If there is in C only one singularity at point a, the value 
of the integral is called the residue of G(s) at the 
singular point a. 

The well known Residues Theorem [Ref. 3] is 

The integral of an analytic function along a 
closed contour is the sum of the residues at 
the singular points in the domain enclosed by 
the contour, multiplied by 2 TT i. 

Now, let us define the analytic function G(s) to be 

G(s) = e . <j6y(s) , 

where (^)-p(s) is the Laplace transform of the sum (convolu- 
tion) T ( + . . . + T n . Then, according to the Theorem of 

Residues, the convolution density will be obtained as 
follows ; 
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f 7 (t) 



= J[0T( s )] 



<* + * 00 



2 JT 



erf- i oo 

orf t i oo 



-y ( s ) ds 



in i 



G ( s ) ds 



o< — / oo 

m 



Uni) f 



Z^ 27r 0 k l 






where m is number of singularities (poles), and Rj^ , 
i=l,...,m, is the residue of the analytic function G(s) at 
the singular point (pole) a ^ . 



Therefore 



f T (0 




5 



that is, the convolution density function is just the sum of 
the residue of the analytic function G(s) at the singular 
points . 

The result established by the last equation would not be 
of much value if it were not possible to compute the residue 
of the analytic function G(s) at the singularities directly 
and without evaluating any integral. Thus, let us now 
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discuss how the residue at a singular point (pole) a of the 
holomorphic function G(s), can be evaluated. 

First the order of the singular point (pole) must be 
established by computing the limit, for s tending to a, of 
(s-a)G(s), (s-a) G(s), and so on, until we find a finite 
limit. The exponent which we find for (s-a) is the order of 
the singular point a. If the singular point(pole) at a is of 
order (multiplicity) n, the residue of the function G(s) at 
a is given by [Ref. 3] 



R 



= lim 
S->- CL 





Therefore, the residue of G(s) at the ( singular point (pole) 
a • is 



R 



= lim 
S -*• Cl i 




d 



ds n '" 





y 



where n j is the order (multiplicity) of the singular point 
(pole) a j , i = 1,..., m. 

By substitution, the convolution density function is 
given by 

m 

f T = Y, R i 
i - 1 
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m 



■E 



i*l 



lim 
S ■+«*,' 



>i-0 1 



^((S-O-O'SCS)) • 



This implies that the convolution density can be obtained by 
limits and differentiations which are easier operations to 
deal with than the integrations. 

Now, we are ready to deal with the three different 
convolution cases which have been mentioned before. we will 
try first to obtain the convolution density for every case. 
Then, the survival functions will be easily obtained by 
integration . 



E. SURVIVAL FUNCTION OF A CONVOLUTION 

In this subsection, we will use the results obtained in 
Subsection D to find the survival function for all the 
possible convolution cases which have been mentioned before. 
This will result in a general equation for the survival 
function of any convolution of independent and exponentially 
distributed random variables. 
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1. A System Having n Dissimilar Failure Rates 



Consider a system defined by the shorthand notation 



where 



EXP ( X , ) + • • • + EXP ( \ n ) 



X; t Xj , i £ j , i,j =1, ... ,n. 



The life time of the system is 



T = T, + ... + T n , 

where Tj is an exponentially distributed random variable 
with a failure rate Xj > i=l, . . . ,n and T| j • • • y T n are 
independent. The Laplace transform of the convolution is 
given by 









) 



If we define the analytic function 

S t , 

G ( s ) = e . (pj (s) , 

then, according to the Theorem of Residues, the convolution 
density is given by 
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f T (t ) = 



n 

i 

1=1 



where n is the number of singularities (number of distinct 
failure rates), and R is the residue of the analytic func- 
tion G(s) at the singular point (pole) -X^ , i = 1,..., n. 

It is clear from the Laplace transform function 
(^)-j-(s) that there are n singularities, each of which has 
an order (multiplicity) one; that is, n | = 1 for i = 1, . . . , 
n. Therefore, the residue of G( s) at the singular 
point (pole) -X; is given by 






R ; = lim 

s— -X;L (n; - 0 i As 



, n ' 1 m; 

Lp((s-<-Xi)) <J(S)) 



= lim 

s 



(s + X; ) G (s ) 



By substitution, the convolution density will become 



f T (t) 



I 



i =1 



R • 



■n 



) lim v (s + 

s— X; 



W'tt, X, 



1=1 



TJ(. A ' ) 
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The survival function of the system will be obtained 
as follows ; 

? T (t) = P ( T > t ) 



OO 




which is the well known formula for the survival function 
for a system having n dissimilar failure rates. 
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2. A System Having n Identical Failure Rates 



Consider a system defined by the shorthand notation 
EXP ( X ) + • • • + EXP ( X ) ( n times ) 

The life time of the system is 

T = T, + ... + T n . 

where Tj , i = 1 n, is an exponentially distributed 

random variable with a constant failure rate X , and T ( , . . . , 
T^ are independent. The Laplace transform of the convolu- 
tion is 

n 

0 T ( S ) = V 

r <=i 




Then 



G(s) 





si-XJ 
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In this case the analytic function G(s) has one sigular 
point (pole) of order (multiplicity) n at the point s = -\ . 

According to the Theorem of Residues, the residue of 
the analytic function G(s) is given by 



n- I 



r = / / m 

S — 



= Urn 
s 



(n- 0! 




1 


n- 1 

d / 


(n - i) ! 


ds n -'V 


i 


rt-/ 

J 






r\ 



C st X 



s + 



xr 



= l ‘ m 

s — ^ - X ^ C 10 — 1 ^ i ds 



c f n i 

— ( e x ) 



X 



n 



(n- i) ! 



n-i -X t 

t e , t>o. 



Therefore, the convolution density is given by 
fy(t) = R 




_Xf 

. e 



- kt 

e , t 0, 
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which is the density function of the Erlang (n , X ) 
distribution. 



Now, the survival function is given by 
F T (t) = P ( T > t ) 



co 



J f (u 



) du 



oo n 



■./ x 



-\u 



(n-l)l 



u 



du 



Integrating by parts yields 





A -Xt 




F T (t) = 


£ (-0! ' 


3 t ^ 0 , 


which is the well 


known formula for the 


survival function of 


the Erlang (n , ^ 


) distribution. 




3 . A Sytem 


Having Some Similar 


and Some Dissimilar 


Failure Rate 




The shorthand notation for the 


system is 


EXP ( X, ) 


+ . . . + EXP ( X, ) 


(n ( times) 


+exp(X 2 ) 


+ . . . + EXP ( X z ) 


( n 2 t ime s ) 
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(n m times ) 



+ EXP(X m ) + • • • + EXP(\ m ) 

where m is the number of distinct failure rates, and n is 
the multiplicity of the i-th failure rate \j , i = l,...,m. 
The Laplace transform of the convolution is 





Xi 

S +• Xi 




Then 



G (s ) 




S t 




In this case, the system has m distinct failure rates, so 
that the analytic function G(s) has m singularities. 
According to the Theorem of Residues the convolution density 
function is given by 

m 

f ( t ) = ^ R j , 

i - I 
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where R is the residue of G(s) at the singular point (pole) 

-X; , i= 1, . . . , 



m . 



Since n j is the multiplicity for the failure rate 

of 

-X; 



Xj, i = l,...,m, it follows that the the residue of the 



analytic function G(s) at the singular point (pole) - /\ | is 
given by 



Rj = lim 



-X ; 









= lim 



i d 



_X.Uni-0! Js 






m 



lim 

S — 



t 

1 


YV'nV^ j Ai 


1 ( n i - 0 ! ds^'A 





rn 



(17 xD 



(n,- - 1 ) ! s-» -Xj 



m 



lim 



i n i“' s t T~T 

S £ + Xj) 

05 1 



n .‘-' */-/ 

Let D = d/ds . Then D = d / ds and 






71 ; 



R 



Urn 



i=l 1 '(",-0! s— \ f r j J' ; 



D 



ri:-l _ 5 £ 



> 7 ? 




- 77 J 
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Let us define the function 



m 



H 



- n. 



i(s) = I I (S * Xj ) - 

j* 1 ' 



Then 



Hi-I 


S f 


S f 


ttj - \ ^ 


D 


e H | ( s ) 


= e 


(D + t) H j ( s ) 

■k - 



and 



rn 



(TT\ ) 

R: = — lim 

- 0 ! S_.Xi 



cf flj-l 

e (D+f) w , (s) 



From the Binomial Theorem the term (D + t) can be repre- 

sented as 
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Let us define the term 




Then 



m n; 



n.-i 






K = 0 



cc h<) + ' - X-,t 

— t e 

«! (rij-K-i)! 



Since 



f T (t) 




5 



it follows that the convolution density is 



m 



m n; -i 



)I I 



1=1 K=° 



CP',*) 

* ! (n;-K-i)j 



-X;t 

t e , t>o. 
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This equation is a general equation for the density function 
of any convolution of exponentially distributed random 
variables . 

Now, the survival function can be obtained as 

follows ; 

F T (t) = p( T > t ) 




(u) 



du 



m ^ _\iU 



{V\-)YY , f u : 

\ i=t / Lj Lj k , 1 J 



J 



u. 



|=l K -0 



"■ - — cO,K) V- (Xit) 



- *' z '' X;t 






,=l k=o 'X 1 7=1 



m r) . m Tij— l n;-K 



2-1 



/TTX "W'V’ V " 1 (Xjfl -Xj t 

= (!=A' )Y j LL ■ - ' e > ^ 0 



/ = ! K-0 2=1 



! XT (*- 0 ! 



This is a general equation for the survival function of any 
convolution of exponentially distributed random variables. 
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In the next subsection, we will give three examples 
of the application of the equation. These examples represent 
all possiple convolution forms. 

F. EXAMPLES ON THE GENERAL EQUATION 

1 . A Convolution of Dissimilar Fai lure Rates 

Using the general equation which was derived in the 
last subsection, this example illustrates how to compute the 
survival function for a system defined by the shorthand 
notation 

EXP ( \| ) + EXP(\ 2 ) h. EXP(\ 3 ) 

where 

X ) t X 2_ > X \ t X3 > X 2. ^ X 3 . 

Since all the failure rates are dissimilar, it follows that 
the total number of distinct failure rates is 

m = 3 . 

The multiplicity of the first distinct failure rate X| is 

n | = I- 

The multiplicity of the second distinct failure rate X^ is 
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The multiplicity of the third distinct failure rate \ is 



n 3 “ 1* 



The survival function of the system is given by 



3 n . 3 "i-l fir* 



z-i 



f t (t, -[nOZV Y e Xit 

I \ (-1 l . — ' £ I \ n i - < ✓ j\| 

i - * k=o z-\ 'A; ( z - u . 



= suml + sum2 + sum3 , 



where suml is the sum at i=l, sum2 is the sum at i=2 and 
sum3 is the sum at i=3. 

At i=l, n, = 1, k has only the value 0, and z has 
only the value 1. The value of the derivative term C(1,0) is 
given by 



C(1,0) 




'( S + X 2 ) ( 5 + \$ ) s = - X , 



(X 2 -X l )(x 3 -X 1 ) 
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Therefore 



suml 



l-l i-K 



^-Z ^3 Z! £ 

K—0 fc-' 



co,<) (X,t) -K* 

*! x;- k c*-o < 



^•1 ^2 ^ 3 



CO, o) (X, 0° e - t 
0! X, 0 ! 



X | X 2 X 



3 




X^ X 3 

(h-KXh-K) 





Similarly, at i=2, =1, k has only the value 0, 

and z has only the value 1. The derivative term C(2,0) is 



and 



sum2 



' ( X,- X 2 )(X 3 - x 2 ) 

x, X 3 e X; 
( X|- X 2 )( Xj- X 2 ) 






Also, at i=3, nj = 1, k has only the value 0, and z 
has only the value 1. The derivative term C(3,0) is given by 



C(3,0) 



( ^-1 ^3 )( ^2 ^3 ) 
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; 



and 



sum3 = 



X. X 



\ /N 2. 



-X 3 t 



( X, - X 3 ) ( X2. - X3) 



By adding suml, sum2 and sum3 , the survival function 
of the system is 



F (t) = suml + sum2 + sum3 



\ ~X 2 t s \ -X 3 t- 

+ X) X 3 e i \ , , e , t >0 



X 2 X 3 e A» A 3 e " . Ai /\2 e. 



(X 2 - Xj)(X 3 -X,') (X,-X^(X 3 -X^) (X,-X3 S )(X 2 -X^ 



2. A Convolution of Identical Failure Rates 

This example illustrates how to compute the survival 
function for a system defined by shorthand notation 



EXP ( X ) + EXP ( X ) + EXP ( X ) • 



This convolution has only one distinct failure rate 
with multiplicity 3. Therefore m=l and n=3. The survival 
function of the system can be obtained as follows; 






/—I K-0 Z = l 



* K! x? (z-0 ! 



, t*0. 
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? y y* cv,k) (Xt) - ^ { 

kb, K! )T 

Since m = 1, and n| = 3, it follows that i = 1 and k 
0,1,2. When k = 0 the derivative term C(1,0) is given by 



C(1,0) 




[ 





L (S+9)' 




X 



= I 




= 1 . 



it follows that the k-th derivative 




0 , k 1. 
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That is, C(l,k) = 0, for k = 1,2 . 
tion for the system will become 



Now, the survival func- 



f t (t) 



X I 



C ( 1 , 0) 

X 3 



(Xtf' e -Xt 

(*-')! 



A -Xt 

A ( Z .,)| e 

z-\ ' 



, t ^ 0, 



which is the survival function equation of the Erlang( 3, A) 
distribution. 



3 . A Convolution of Similar and Dissimilar Failure 
Rates 

This example illustrates how to compute the survival 
function for a system having some similar and some dissim- 
ilar failure rates, using the general equation. 

Suppose we have a system defined by the shorthand 

notation 

EXP ( ) + EXP(X j ) + EXP ( ) 

In this case, we have two distinct failure rates X| and 
\ . The total number of the distinct failure rates is 

m = 2 
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The first failure rate has the multiplicity 



n j — 2. « 

The second failure rate X^ has the multiplicity 

n 2 _ = 1 • 

The survival function equation of the convolution can be 
given as follows; 



2 Hr- 1 *;-< 



, z-i \ 

c ci, <) (X; t) -X/t 

„/ \ r*'" K rr. O' 

i=l Z=M • A.j V. Z- '/ . 






, t > 0 



= suml + sum2 , 

where suml is the sum at i=l, and sum£ is the sum at i=2. 

At i=l, k takes the values 0 and 1. We will evaluate 
the two derivatives C(1,0) and C(l,l) as follows; 



C(1,0) = 



fr U*^jY J 



s jt ' 

2 



.-X, 



U( s+ h) J 



j 



s = -X, 
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U 2 -X,) 



C(l,l) = 



d s' i 



xix * 

J? 1 






-X, 



(x 2 -x,) : 



Then 

suml 



2 -< 



2 — > I 



-X, t 



■ x* x z £ I e 

2 fe»z=, «! X 2 (2-1) ! 

Z , C<l,o) rXi't CO,o)X L t -X,t CO,0 r X|t ' 



= X, X 



I A 2 



x; 



e ± 



X; 



e 1 1 



X, 



Xz X 1 Xzt X, Xz 



Xz \ x 2 - x, (x 2 -x,y 



x,t 



At i = 2, n^ =1 and k takes only the value 0 . 
we have to evaluate only the derivative term C(2,0) 
follows ; 





r d° 


TT/ x i 


C ( 2 , 0 ) = 


i o 


WO + Xj) 




d s 


J>2 J 



-X: 



(X r X 2 ) ; 
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Thus , 
as 



Then 



2 . 

sum2 = X, X z 



C(2,°) 

0 ! \ 2_ o ) 



k ' Xi (X t -X ! ) 1 X ! * 

^ e -M . 

( x*-x,) 



Now, the survival function for the system is given by 



F 



T 



(t) 



suml + sum2 



Xit 

- Xz~ X| X; 




The three examples above illustrate how we can 
compute the reliability of any convolution of exponenial 
random variables. In practice, a computer program computes 
the reliability directly. A practitioner of the shorthand 
methodology can use the computer program in Appendix A to 
obtain a system's survival function by just inserting the 
failure rates of the system components. 
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III. SURVIVAL FUNCTION AND POLYNOMIAL COEFFICIENTS 



In section 2, we derived a general equation for the 
survival function of any convolution of independent and 
exponentially distributed random variables. This section 
gives an alternative formula for the survival function, 
using the idea of computing the coefficients of the polyno- 
mials that accompany the exponential terms in the survival 
function equation. This idea was explored by Sadan Gursel 
[Ref. 4], adding one distribution at a time. In this section 
we will derive the same polynomial coefficients from the 
general equation developed in the previous section. 

A. DETERMINATION OF THE POLYNOMIAL COEFFICIENTS 

Consider a system defined by the shorthand notation 
EXP ( \ j ) + ... + EXP ( X | ) (n, times) 

+ EXP (X ^ ) + • • • + EXP(\ Z ) ( n 2 times) 

+ 

+ EXP(\^ )+ ... + EXP ( \ ^ ) (n^ times). 
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m 



This notation represents the convolution 
tial random variables, where there are n 
tial random variables having the failure 



I nj 



of 

i r » 

identical 



rate 



exponen- 

exponen- 



X | , i = 1 , ••• , m . 

We expect the survival function of the convolution to 
have the form 

_X,t -X 2 t -X m t 

F (t) = A ^ ( t ) . e + A^(t ) . e + ... + A m (t).e , t £ 0, 



where 



A j ( t ) = 



10 



\\ 



IZ 



2 

t + 



n r i 



i n.- \ 



z 

A 2(t) = a 2 o + a 2i + a 22 t + 



+ 



a 2n r l 



t 






2 n rr|" 1 

A m( t ) = a mo + a mi t + a mz t + • • • + a nvn rT p ' t 

The notation Aj(t), i=l,...,m, represents a polynomial of 
(nj-l)st degree. Therefore, to obtain the survival function 
of the convolution, we have to find all the coefficients 
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For example, 
notation 



the convolution represented by the shorthand 



EXP(\,) + EXP ( \ | ) + EXP ( ) + EXP ( \ 2 ) , 

where 

X | / ’ 

has a survival function represented by 



F (t) = A , ( t ) 




+ A 2> (t) 




t > 0, 



where 

Z 

A | ( t ) = a ( g + a,j t + a t , 

and 

A 2 (t) = a £o • 

It is important to note that the number of exponential 
terms is equal to the number of dissimilar failure rates, 
and each exponential term has a polynomial coefficients of 
degree one less than the number of identical random vari- 
ables having the corresponding failure rate. 

The problem of computing the coefficients will be solved 
if we derive a general formula for the constant coefficients 
of the polynomial coefficients in the survival function 
equation. 
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Let us begin with the general equation which is 



m y,. m rtf- 1 n;-K v . . 



-X;t 



, t£0 



)rt 



rt; 



The term TT \! can be written as 
hi I 



TT h 1 ^ n J 

w X| ■ x i(uV) • 

J* > 

Substituting the above into the general equation yields 



Z-l 



' <■> -EC( f J Xi ' 

i = f jti K*o Z*\ h[ (Z-ni 



Rearranging the above equation yields 



X m n;-i . tr+z-/ z-l 

^(ITXj )y y £<i^Lx. t 



Let L = z - 1. Then the survival function will become 



m \ + 70 n; n r i nri-K k + l I 

'( 1 ME frir^i ( 

i -i J>» k=oL= ; o 
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Let us define 



K + L 

G,(k,L) = COjjO X t 

K > LI 



Then, 






Tl ( -) (ft;- \)-< 


L L G 


= G; (0,0)+ 


+ G j (0 , n 


K-0 L = 0 


+ G j (1,0) + 


+ G j ( 1 , n- - 2 ) 




+ G j (2 , 0 ) + . 


• - + Gj (2,n; -3) 



+ 



+ 

+ ... 

+ G | (n r l,0) 

I 

or, equivalently 



^ G j ( k , L ) = G j (0,0)+ +G;(n r l,0) 

K =° L = C> + G j ( 0 , 1 ) + + G (n j -2,1) 

+ G-, (0,2)+. . . +G j ( n • - 3,2) 



+ G j ( 0 , n- - 1 ) 
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rij-i (n ; -i)-L 

Y Gj(k,L). 

1=0 K = 0 

Therefore 




K+- L L 



m -\.t m t\- *i-l . v 

T «"I« ' (TXj J >I I t 

JVi “> «■• ' l - 



m 



r~> — \ . "t 

^ e A l . A . ( t ) , 



i=l 



where 



*;-/ 



V‘> ■ Z (TT,^)! 






L = 0 ff; K=0 



CU 2 *) 

K ! L! 



k+ L L 

X; t • 



Now, it easily follows that 



a . 



m i\\ fa-d - L \ K+ h 

C Q, *0 



I L 



- ( I 



J sl 



K=0 



X, 

K! LI ' 



represents the coefficient of t for the i-th polynomial in 
the survival function equation, L = 0, ... , nj-1, 
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where 




This formula is the general formula fqr the coefficients and 
it can be considered as an alternative formula for the 
survival function of any convolution of exponential random 
variab les . 

It is clear from the above formula and from the general 
equation in Section 2 that the most difficult part in the 
reliability computations is the derivative term C(i,k), 
particularly when dealing with complex convolutions that 
have a long stream of failure rates with high corresponding 
multiplicities. We will give the the practical method for 
computing the term C(i,k) in Subsection C. 

B . EXAMPLE 

Consider a system defined by the shorthand notation 
EXP ( \ , ) + EXP ( X ^ ) + Exp (\ 3 ) + EX p (X 3 ) 

Then, the survival function for the system is given by 

F (t) = A, (t).e“^ |1: + A a (t).e"^ t + A 5 ( t ) . e 



51 



The first distinct failure rate \| has a multiplicity 1, 
the second distinct failure rate has a multiplicity 1, 
and the third distict failure rate X3 has a multiplicity 2. 
Therefore 



m = 3 , 
n, = 1 , 
n 2 = 1 , 

n 3 = 2 . 

The polynomial A|(t) will be obtained from 

" 1 - 1 

A, (t) = a [0 + a n t + . . . + a |r1r , t 
Substituting nj = 1 in the above equation yields 
A ,(t) = a )0 . 

Since 



1 L 



m n- 

( TUi ) I 



rt ; -1 -L 



j-'. 

J & i 



k=o 



c c i,<) 

Kt LI 



K + L 

X,- , 



it follows that 



a l 0 






l-l-O 



J=l 



K=o 



CO, *0 
<! 0 ! 



K + 0 

X, 
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c(/> *0 

K 1 



Z 0 

- x 2 x 3 y . 

K = 0 




z 

=X 2 X 3 c d.°) . 



where 



C(1,0) 







l 

(Xz-X,)(X 3 -X l ) 2 



Thus 



2 

Xz X3 

,0 = (X 2 -X,)(X 3 -X 1 ) Z 



Now, the polynomial A|(t) is 
A | (t) = a | 0 



Z 

Xji X5 

(X 2 -X,)(X 3 -X ,) 2 



Similarly, the polynomial A^(t) is given by 
A 2.( ** ) = a 2. 0 
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X, X 



t A 3 



(X, X 2 ) ( X 3 X 2 .) ' 



The polynomial A 3 (t) will be obtained from 



A 3 (t ) = a 3Q + t 



+ a 3 rij - 1 t 



n 3“ l 



Substituting n 3 = 2 yields 



A 3^) = a 30 + a 3l * 



In this case we have to evaluate the two coefficients 
a^ and a 3 ( . Substituting nj =2, m = 3, L = 0 in the 
general equation of the coefficients yields 



30 



■(^\i )E 

J #3 K = 



■ -"° C ( 3 , k) 



K = 0 
1 



x, X2 Y 



K { 0 ! 

C (3 , K) 



K = 0 



K ! 



x. 



X, 



K +0 



= X| X Z C( 3 , 0 ) + X | X^C ( 3 , 1 ) X 3 , 



where 



C ( 3 , 0 ) = 



ds° 



3 - n j 

n<> + Xj) 

J>3 



S = -X 5 



(X,-X 3 )(X 2 -X 3 ) ’ 
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and 



_ r ) 



C ( 3 , 1 ) = 



Aft* + ) 

dS J>3 



J 1 



-X 3 

I 



' ( Aj Aj) ( X^- Xj) (Aj~ Aj)(X2“\s) 



Thus 



X, X 



I ^2 



A 1 A? A 



2^3 



A ) A 2 A 



1 ^2^3 



*30 



(ArAsXAz-)^) (A, A 3 )( A2~ A3) (Aj A 3 )(A2 A 3 ) 



of coefficients yields 



\31 



j *3 



X,Xz t 



: 2 , and 


m = : 


2-1-1 


C(3 
K J 


) I 


K=0 




C (3 


7 ^ 



< = 0 



K ! 



A, A? c ( 3 , 0 ) x*-z . 



Since 



3 

3 



C (3 , 0 ) = 



K+* I 
•3 



X: 



< + * 



r 1 — r _7 V 

[ 77 j) 

J>3 



S ’" " A 3 
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( X,- X 3 )( Xr X 3 ) 



it follows that 



a 3 / 




(X.-^XXz-X,) 



Finally, the survival function of the convolution is 
given by 



F T (t) =A,(t).e^'^ - A A A 3 



+ A^(t ) .e + A j ( t ) e ^ 



-X,t 



-X?t 



= a )0 e A l + a 7n e 2 + ( 



2.0 



a 30 +a 3l 



-X,t 

t ) e -5 



x.xU^ 



2- - X Z t 
X| X3 e 



( X 2 ~ X,)( X 3 - X,) (X, X 2 )(X 3 X 2 ) 



Xi X 2 



X, X 2 X3 



KXf X 3 )( Xr X 3 ) (x,-xj(x 2 -x 3 ) 



X, Xz Xj 



( X| X 3 )( Xz 



X, X 2 X 3 1 j ^-Xjt 



r X 3 )(X 2 X 5 ) 
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C. THE PRACTICAL METHOD FOR COMPUTING THE DERIVATIVE TERM 



C(I,K) 

One of the major difficulties in computing the survival 
function either by the general equation, or by the general 
formula of the coefficients, is the derivative term C(i,k), 



where 



C ( i , k) 





m 




j* * 








) 



This term represents the k-th derivative for a product of 
m-1 terms. The difficulty in comuting C(i,k) increases with 
k. Since k runs from zero up to nj -1, i=l,...,m, it follows 
that the maximum value of k will not exceed the maximum 

I 

multiplicity of the distinct failure rates minus one. For 
example, if we have the two distinct failure rates \ and 
X z where has a multiplicity n^ and has a multiplicity 
n^ , n ( > n^ , then, the maximum value of k is n ( - 1. 

In practice, due to the high accuracy (precision) 
required in the computations, we will limit the derivative 
term C(i,k) to values of k <9. This implies that the 
multiplicity of any of the distinct failure rates in any 
convolution must not exeed 10. 
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The following is one of the solutions which will be used 
in the computer program for computing the term C(i,k) (see 
Appendix. A, Appendix. B) . 

Suppose we have a convolution which has m distinct 
failure rates with in corresponding multiplicities . We know 
that 



Suppose for example , we want to compute C(i,k) at i = 1 , 
then k=0, ... , n ( -l, and 



Our problem is how to change a derivative of a product of 
terms to a derivative of a sum of terms, since the later is 
easier in computations. 

Let us define 




m 




m 



C(l,k) 





m 



f 



j*' 
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At i=l we obtain 



ffl 



- n 



= II ( s ^j) 

j * 1 

i 

= (T7xJV7 Jj r ~(77x^r 



By taking the natural logarithiin on both sides, we obtain 



In f | = - 



rt 2 



ln(s + + ••• + ln ( s + X,„) 



>V 



Now, by taking the first derivative on both sides , 
obtain 



f, / f , = 



^ 2 . 



5 + X; 



rb 



^ + Xm 



or, equivalently 



y 
f , 



h 



s + X; 



+ 



r >, 



^ Xyn - 



Let us define the term 



m 



-I 



-n 



J = ' 

j * ' 



\j 



we 
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Then , 



m 



A i = 



I 

j=' 



n J 



+ Xj 



r\ z 



n 



m 



S + 



X : 



-h 



S t 



X 



m 



Then 



i 




or, equivalently 



df j /ds = - A | . f j . 

When i=l and k= 1 , the term C(i,k) is 



C(l,l) 




m - rij 

ljJ( £ + \i ) 



<J 5 




But 
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so that 



C(l,l) 




We can proceed in a similar manner to compute the 
derivatives C ( 1 , 2 ) , C ( 1 , 3 ) , . . . , C ( 1 , n -1) as follows; 



where 



f , = df, / ds 



d/ds 



- A, . f. 



but 



so that 



- - A, . f, - A, . , 



f , - - Aj . f| , 



f, 



= ‘ A l f l + A | f l 



zr z 



= (-D 



A | " A, 



. f 



rr\ 



A, 



= d/ds ^ 



n 



j* 1 



s-t- Xj 



other 
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TT) 



I 



j*' 

jf' 




Then, the derivative term C(l,2) is 
C ( 1 j 2 ) = 



~X, 



, .2 


r 2- / i 


(-0 


( A| - A, } f, 
m * 



--X, 



and so on until we obtain the derivative term C(l,nj-1) (see 
Appendix . A , Appendix . B ) . 
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APPENDIX A 



This appendix contains a computer program written in 
Fortran for the survival function of any convolution of 
independent and exponentially distributed random variables, 
using the general equation 



"■>-(m 1 ')EE f^E 

i=l • i is-n ! A ; _ 



ir,-K 



z-l 



i=i K =o 



K! A, 



z-i 



(Ajt) 

(Z-O ) 



Ait 

e , t >o. 



A. FUNCTIONS AND SUBROUTINES USED IN THE PROGRAM 

1. FUNCTION I FAC 

This function computes any factorial required in the 

program. 

2. SUBROUTINE GAMMA 

This subprogram computes the incomplete gamma func- 
tion with the parameters( ,nj-k) and returns the result in 
the variable G, 
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where 




3. SUBROUTINE DERIV 

This subprogram computes the k-th derivative C(i,k) 
for k<9 and returns the result in the variable C. 

B. PROGRAM'S LIMITATIONS 

(1) The multiplicity of any of the distinct failure rates 
must not exceed 10. 

(2) The computer precision limit. 

(3) There is no limitation on the number of distinct failure 
rates . 
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nonnnnooooooo 



section Computer program 



$ JOB 

LAMDAI . . . . I - th distinct failure rate. 

LAMDA Array of M distinct failure rates. 

PRODCT .... Product of all the failure rates. 

FBAR Survival function of the convolution. 

T Time. 

C Value of the derivative term C(i,k). 

G Incomplete gamma function. 

MULTI Multiplicity of the i-th distinct 

failure rate. 

MULTIP ... .Array of multiplicities of the distinct 
failure rates. 

M # of distinct failure rates. 

I,J,K Loop index. 



DOUBLE PRECISION LAMDA (20), LAMDAI , PRODCT , FBAR , 

*T , C , G 

INTEGER I , J ,N ,M ,K , MULTI , MULTIP (20 ) 

PRINT, 'Please enter the time of duration' 

READ, T 

PRINT, 'Please enter THE # OF dissimilar failure rates' 
READ, M 

PRINT, 'Enter the values of dissimilar failure rates' 
READ, ( LAMDA ( J ) , J=1,M) 

PRINT, 'Enter corresponding multiplicities ' 

READ, (MULTIP ( J ) ,J=1,M) 

PRINT , ' OUTPUT ' 

PRINT, '======' 

WRITE ( 6 , 11 )T 

11 FORMAT (' TIME = ’,F10.5) 

PRINT,' FAILURE RATE MULTIPLICITY' 

FBAR = 0.0D0 
PRODCT = 1.0D0 
DO 10 I = 1 ,M 

WRITE (6 , 12 ) LAMDA (I ) , MULTIP (I) 

12 FORMAT (8X,F10.5,10X,I3) 

MULTI = MULTIP (I) 

PRODCT =PRODCT*LAMDA( I ) —MULTI 
DO 20 J=l, MULTI 
K = J - 1 
LAMDAI = LAMDA ( I ) 

CALL GAMMA ( MULTI, K,T, LAMDAI, G) 

CALL DERIV( I, M, LAMDA, MULTIP, K,C) 

FBAR=FBAR+C*G/ (IFAC(K)*LAMDAI** (MULTI-K) ) 

20 CONTINUE 
10 CONTINUE 

FBAR = FBAR * PRODCT 
WRITE (6, 13 )T, FBAR 

13 FORMAT ( IX, 'P(T >',F6.2,') = ',F10.5) 

STOP 
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30 

40 

50 

70 

60 



END 

FUNCTION IFAC(N) 

JFAC= 1 

IF(N.NE.O) THEN DO 

DO 30 J= 1 , N 

JFAC=JFAC*J 

CONTINUE 

END IF 

IFAC= JFAC 

RETURN 

END 

SUBROUTINE GAMMA (MULTI , K , T , LAMDAI , G ) 

DOUBLE PRECISION T , G , LAMDAI 
INTEGER L,K,Z, MULTI 
G= .0D0 
L=MULTI-K 
DO 40 Z=1 ,L 

G = G+ (LAMDAI *T ) ** ( Z - 1 ) *DEXP ( - LAMDAI *T ) / IFAC ( Z - 1 ) 

CONTINUE 

RETURN 

END 

SUBROUTINE DERIV (I ,M , LAMDA , MULTIP ,K , C ) 

DOUBLE PRECISION A , F , C , B ( 9 ) , LAMDA ( 10 ) 

INTEGER I , J ,M , L , IK , MULTIP ( 10 ) 

A= . 0D0 
F= 1 . DO 
DO 50 J= 1 , M 
IF(I.NE.J) THEN DO 

F = F* ( LAMDA ( J ) - LAMDA ( I ) ) ** ( - MULTIP ( J ) ) 

A= A+ ( MULT IP ( J )/( LAMDA ( J )- LAMDA ( I )) ) 

END IF 
CONTINUE 

IF (K . GT . 1 ) THEN DO 
IK=K- 1 

DO 60 L= 1 , IK 

B (L ) = . 0D0 

DO 70 J=1,M 

IF(I.NE.J) THEN DO 

B ( L ) = B (L )+(((- 1 ) **L ) '-MULTIP ( J ) * 

* ( ( LAMDA ( J ) - LAMDA ( I ) ) ** ( - L - 1 ) ) )*IFAC(L) 

END IF 
CONTINUE 
CONTINUE 
END IF 

IF(K.EQ.O)THEN DO 
C = F 

ELSE DO 

IF (K . EQ. 1)THEN DO 
C= (-1)- V A“F 
ELSE DO 

IF (K. EQ. 2 ) THEN DO 
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C=((A**2)-B(1))*F 
ELSE DO 

IF (K . EQ . 3 )THEN DO 

C= ( ( - 1 ) * (A**3 ) + 3*A*B ( 1 ) -B ( 2 ) ) *F 

ELSE DO 

IF(K.EQ.4)THEN DO 
C = ( A**4 - 6 *a**2*B ( 1 ) + 4 * A*B ( 2 ) + 3 * 

*B(1)**2-B(3) )*F 
ELSE DO 

IF(K.EQ. 5 )THEN DO 

C= ( ( - 1 ) * (A** 5 ) + 10* (A**3 ) *B ( 1 ) - 10* ( A**2 ) * 

*b ( 2 ) - 15 *A*B ( 1 ) **2 + 5 *A*B ( 3 ) + 10*B ( 1 ) *B ( 2 ) - B ( 4 ) ) *F 



ELSE DO 

IF (K . EQ. 6 )THEN DO 
C = ( A** 6 - 15 * ( A**4 ) *B ( 1 ) + 2 0 * ( A** 3 ) * 

*B ( 2 ) - 15 * ( A** 2 ) *B ( 3 ) + 6 *A*B ( 4 ) - B ( 5 ) + 4 5 * 

* ( A**2 ) * (B ( 1 ) **2 ) - 6 0 . DO*A*B ( 1 ) * 

*B ( 2 ) + 1 5 *B ( 1 ) *B ( 3 ) - 1 5 *B ( 1 ) ** 3 + 10 *B ( 2 ) **2 ) *F 



ELSE DO 

IF (K . EQ. 7 )THEN DO 

C= ( ( - 1 ) * ( A** 7 ) + 2 1*a** 5*B ( 1 ) - 35*A**4*B ( 2 ) + 35*A**3* 
b ( 3 ) - 21*A**2*B (4 ) + 7*A*B ( 5 ) -B ( 6 ) - 105*A**3* 
b ( 1 ) **2+ 2 10*A**2*B ( 1 ) *B ( 2 ) - 105*A*B ( 1 ) *B ( 3 ) - 105* 

B ( 1 ) ** 2 *B ( 2 ) - 7 0 * a* 

B ( 2 ) **2+ 10 5*A*B ( 1 ) **3 + 2 1*B ( 1 ) *B ( 4 ) + 3 5 *B ( 2 ) * 

B ( 3 ) ) *F 
ELSE DO 

IF (K . EQ. 8 )THEN DO 

C= (A**8 - 28*A**6*B ( 1 ) + 5 6*A**5*B ( 2 ) - 70*A**4*B ( 3 ) + 5 6* 
*A**3*B ( 4 ) - 2 8 *A**2*B ( 5 ) + 8 *A*B (6)-B(7)+210* 

*A**4*B ( 1 ) **2 - 5 60*a**3*b ( 1 ) * 

*B ( 2 ) + 4 2 0 *A** 2*B( 1)*B (3 ) - 168 * A*B ( 1 ) *B ( 4 ) + 2 8 * 

*B ( 1 ) *B ( 5 ) - 420*A**2*B ( 1 ) **3 + 840*A*B ( 1 ) **2*B ( 2 ) + 280* 
*A**2*B ( 2 ) **2 - 2 10*B ( 1 ) **2*B ( 3 ) - 28 0*A*B ( 2 ) *B ( 3 ) - 28 0*B ( 1 ) 
*B ( 2 ) **2 + 5 6 *B ( 2 ) *B ( 4 ) + 3 5 *B ( 3 ) **2 + 1 0 5 *B ( 1 ) **4 ) * F 
ELSE DO 

IF (K . EQ . 9 )THEN DO 

C = ( - A** 9 + 36 *A** 7 *B ( 1 ) - 8 4*A** 6 *B ( 2 ) + 12 6 *A**5 *B ( 3 ) - 1 2 6 * 
*A**4*B ( 4 ) + 84*A**3*B ( 5 ) - 3 6*A**2*B ( 6 ) + 9*A*B ( 7 ) -B ( 8 ) - 
*378 *A** 5 *B(1) **2+1230 *A**4*B ( 1 ) *B ( 2 ) - 1 2 6 0 *A** 3 *B ( 1 ) * 

*b ( 3 ) + 7 5 6 *A**2*B ( 1 ) *B ( 4 ) - 2 5 2*A*B ( 1 ) *B ( 5 ) + 3 6 *B ( 1 ) *B ( 6 ) ) 
C=C+ 1260*A**3*B ( 1)**3- 3780*A**2*B ( 1)**2* 

*B ( 2 ) - 840*A**3*b ( 2 ) ** 2 + 18 90*a* 

*B ( 1 ) **2*B ( 3 ) + 12 6 0*A**2*B ( 2 ) *B ( 3 ) - 3 7 8*B ( 4 ) * 

*B ( 1 ) **2 - 404*A*B ( 2 ) *B ( 4 ) + 84*B ( 2 ) *B ( 5 ) - 945*a* 

*b ( 1 ) **4 + 2520 *A*B (1) *B (2) **2+1260 *B(1) ** 3 *B ( 2 ) - 3 1 5 *A* 

*b ( 3 ) **2- 12 6 0*B ( 1 ) *B ( 2 ) *B ( 3 ) + 12 6*B ( 3 ) *B ( 4 ) - 28 0*B ( 2 ) **3 
C = C*F 



END IF 
END IF 
END IF 
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END IF 
END IF 
END IF 
END IF 
END IF 
END IF 
END IF 
RETURN 
END 

$ENTRY 
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section Examples subsection All failure rates are dissimilar 
Shorthand notation : 

EXP ( . 3 ) +EXP ( . 32)+ EXP ( . 4 ) +EXP ( . 6 ) +EXP ( .62) 



Output 




TIME = 


2.00000 

FAILURE RATE MULTIPLICITY 

0.30000 1 

0.32000 1 

0.40000 1 

0.60000 1 

0.62000 1 


A 

H 


2.00) = 0.99817 
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subsection All failure rates are identical 



Shorthand notation : 



EXP (2.1) + ... + EXP (2.1) (10 times) 



Output 



TIME = 3.20 

FAILURE RATE MULTIPLICITY 

2.10000 10 

P (T > 3.20) = 0.85772 
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subsection Some similar and some dissimilar failure rates 



Shorthand notation : 

+ EXP (6.5)+ ... + EXP (6. 5) (5 times) 



EXP (4. 5)+ 


. . . + EXP (4.5) 


(5 times) 


EXP (3.5)+ 


. . . + EXP (3.5) 


(5 times) 


EXP (2.5)+ 


. . . + EXP (2.5) 


(5 times) 


EXP (1.5)+ 


. . . + EXP (1.5) 


(5 times) 


EXP ( .50)+ 


. . . + EXP ( .50) 


(5 times) 



Output 



TIME = 8.00000 



FAILURE RATE 

6.50000 

4.50000 

3.50000 

2.50000 

1.50000 
0.50000 



MULTIPLICITY 

5 

5 

5 

5 

5 

5 



P (T > 8.00) = 0.99894 

cms cms 
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APPENDIX B 



This appendix contains an alternative computer program 
written in Fortran for computing the survival function of 
any convolution of independent and exponentially distributed 
random variables, using the general equation for the coeffi- 
cients of the polynomials that accompany the exponential 
terms in the survival function equation 



F ( t ) = 



I 




L= o 



L 

O-iL t 



t >0 , 



where 



a 



iL 



m 



"j 



( TTTXj ) 

N j=i j ' 
j^i 



I 



K -0 



. K+ L 

Cih <) 

K] LI ' 



A. FUNCTIONS AND SUBROUTINES USED IN THE PROGRAM 
1. Function IFAC 



This function computes any factorial required in the 

program. 
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2. Subroutine COEFF 



This subprogram computes the coeffient of t for the 



i-th polynomial in the survival function equation and 
returns the result in the variable COj 

where 



3. SUBROUTINE DERIV 

This subprogram computes the k-th derivative C(i,k) 
for k<9 and returns the result in the variable C. 

B. PROGRAM'S LIMITATIONS : 

(1) The multiplicity for any of the distinct failure rates 
must not exceed 10. 

(2) The computations must be within the computer precision 
limit . 

(3) No limitation on the number of distinct failure rates. 




J 
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oooooooooooono 



section Computer program 



$ JOB 

LAMDAI . . . . i- th distinct failure rate. 



LAMDA Array of M distinct failure rates. 

MULTI Multiplicity of the i-th distinct 

failure rate. 



MULTIP .... Array of multiplicities of the distinct 
failure rates. 

PRODCT .... Product of the distinct failure rates. 



FBAR Survival function of the convolution. 

T Time . 

C Value of the derivative term C(i,k). 

CO Coeffetiont of the polynomials. 

M # of distinct failure rates. 

I,J,K Loop index. 



DOUBLE PRECISION LAMDA ( 20 ) , LAMDAI , PRODCT , FBAR , 

*T , C , CO 

INTEGER I , J , N ,M , K , MULTI , MULTIP (20 ) 

PRINT, 'Please enter the time of duration' 

READ, T 

PRINT, 'Please Enter # of dissimilar failure rates' 
READ, M 

PRINT, 'Enter the values of dissimilar failure rates' 
READ, ( LAMDA ( J ) , J=1,M) 

PRINT, 'Enter corresponding multiplicities ' 

READ, (MULTIP ( J) ,J=1,M) 

PRINT, 'OUTPUT ' 

PRINT, ' ====== ' 

WRITE ( 6 , 11 )T 

11 FORMAT ( ' TIME = ’,F10.5) 

PRINT,' FAILURE RATE MULTIPLICITY' 

FBAR = 0.0D0 
DO 10 I = 1,M 

10 WRITE (6,12 )LAMDA ( I ) , MULTIP ( I ) 

12 FORMAT (8X,F10.5,10X,I3) 

DO 20 I = 1 , M 

MULTI = MULTIP (I ) 

WRITE (6,13)1 

13 F0RMAT( IX, 'POLYNOMIALS COEFFICIENT FOR LAMDA', 12,' is 
DO 30 LL = 1, MULTI 

L = LL - 1 

CALL COEFF ( I, L,M, MULTIP, LAMDA, CO) 

WRITE ( 6 , 14)1, L, CO 

14 F0RMAT(3X, 'a' ,11,11, ' = ’.F25.12) 

FBAR= FBAR+ T*' V L " CO "DEXP ( - LAMDA ( I ) *T ) 

30 CONTINUE 
20 CONTINUE 

WRITE (6 ,15 )T, FBAR 

15 FORMAT ( / , IX , ' P (T >',F6.2,') = ’,F20.12) 
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40 

50 

60 

70 

90 



STOP 

END 

FUNCTION IFAC(N) 

JFAC= 1 

IF(N.NE.O) THEN DO 

DO 40 J=1,N 

JFAC=JFAC*J 

CONTINUE 

END IF 

IFAC= JFAC 

RETURN 

END 

SUBROUTINE COEFF (I ,L,M, MULTIP , LAMDA , CO ) 

DOUBLE PRECISION C , CO , LAMDA ( 10 ) , PRODCT 
INTEGER I , L , M , K , KK , LL , MULTIP ( 10 ) 

PRODCT = 1.0D0 

DO 50 J = 1 ,M 

IF(I.NE.J) THEN DO 

PRODCT = PRODCT* LAMDA ( J ) **MULTIP ( J ) 

END IF 

CONTINUE 

CO = 0.0D0 

LL = MULTIP (I) - L 

DO 60 KK = 1 , LL 

K = KK - 1 

CALL DERI V ( I , M , LAMDA , MULTIP , K , C ) 

CO = CO + C*PRODCT*LAMDA ( I ) **(K+L ) / IFAC (K) / IFAC (L) 

CONTINUE 

RETURN 

END 

SUBROUTINE DERIV (I ,M, LAMDA , MULTIP ,K , C ) 

DOUBLE PRECISION A , F , C , B ( 9 ) , LAMDA ( 10 ) 

INTEGER I , J , M , L , IK , MULTIP ( 10 ) 

A= .0D0 
F= 1 . DO 
DO 70 J= 1 ,M 
IF(I.NE.J) THEN DO 

F=F*(LAMDA( J ) -LAMDA(I ) )**( -MULTIP (J ) ) 

A=A+ (MULTIP (J) / (LAMDA(J)-LAMDA(I) ) ) 

END IF 
CONTINUE 

IF (K . GT . 1) THEN DO 
IK=K- 1 

DO 80 L= 1 , IK 

B (L) = . 0D0 

DO 90 J= 1 ,M 

IF(I.NE.J) THEN DO 

B(L)=B(L)+ ( ( ( - 1 ) **L ) *MULTIP ( J ) * 

*( (LAMDA(J)-LAMDA(I) )**(-L-l) ) )*IFAC(L) 

END IF 
CONTINUE 
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80 



CONTINUE 
END IF 

IF (K . EQ . 0 ) THEN DO 
C = F 

ELSE DO 

IF (K . EQ. 1)THEN DO 
C= (- 1)*A*F 
ELSE DO 

IF (K . EQ . 2 ) THEN DO 
C=((A**2)-B(1))*F 
ELSE DO 

IF (K . EQ. 3 )THEN DO 

C= ( (-1)*(A**3)+3*A*B(1)-B(2) )*F 

ELSE DO 

IF(K.EQ.4)THEN DO 
C= ( A**4 - 6 *A** 2 *B ( 1 ) +4*A*B ( 2 ) + 3* 
*B(1)**2-B(3) )*F 



ELSE DO 

IF (K . EQ . 5 )THEN DO 

C= ( ( - 1 ) * ( A**5 ) + 10* (A** 3 ) *B ( 1 ) - 10* ( A**2 ) * 

*B ( 2 ) - 15 *A*B ( 1 ) **2 + 5 *A*B ( 3 ) + 10*B ( 1 ) *B ( 2 ) - B ( 4 ) ) *F 

ELSE DO 

IF(K.EQ. 6) THEN DO 

C = ( A** 6 - 1 5 * ( A* *4 ) *B ( 1 ) + 2 0 * ( A** 3 ) * 

*B ( 2 ) - 15* ( A**2 ) *B ( 3 ) + 6*A*B (4)-B(5)+45* 

* ( A**2 ) * ( B ( 1 ) **2 ) - 6 0 . DO*A*B ( 1 ) * 

*B ( 2 ) + 15*B ( 1 ) *B ( 3 ) - 15*B ( 1 ) **3 + 10*B ( 2 ) **2 ) *F 

ELSE DO 

IF (K . EQ. 7 )THEN DO 

C= ( ( - 1 ) * (A** 7 ) + 2 1*A**5*B ( 1 ) - 3 5 *A**4*B ( 2 ) + 3 5 *A**3* 

*B ( 3 ) - 2 1*A**2*B ( 4 ) + 7*A*B ( 5 ) - B ( 6 ) - 10 5*A**3* 

*B ( 1 ) **2 + 2 10*A**2*B ( 1 ) *B ( 2 ) - 105*A*B ( 1 ) *B ( 3 ) - 105* 

*B ( 1 ) ** 2 *B ( 2 ) - 7 0 *A* 

*B ( 2 ) **2 + 10 5*A*B ( 1 ) ** 3 + 2 1*B ( 1 ) *B ( 4 ) + 3 5*B ( 2 ) * 

*B ( 3 ) ) *F 

ELSE DO 

IF (K . EQ. 8 )THEN DO 

C = ( A** 8-28 *A** 6 *B ( 1 ) + 5 6 *A**5 *B ( 2 ) - 7 0 *A**4*B ( 3 ) + 5 6 * 
*A**3*B ( 4 ) - 28*A**2*B ( 5 ) + 8*A*B (6)-B(7)+210* 

*A**4*B ( 1 ) **2 - 5 6 0*A**3 *B ( 1 ) * 

*B ( 2 ) + 420 *A** 2 *B ( 1 ) *B ( 3 ) - 1 6 8 *A*B ( 1 ) *B ( 4 ) + 2 8 * 

*B ( 1 ) *B ( 5 ) - 420*A**2*B ( 1 ) **3 + 840*A*B ( 1 ) **2*B ( 2 ) + 280* 
*A**2*B ( 2 ) **2 - 2 10*B ( 1 ) **2*B ( 3 ) - 2 8 0*A*B ( 2 ) *B ( 3 ) - 28 0*B ( 1 ) 
*B ( 2 ) ** 2 + 5 6 *B ( 2 ) *B ( 4 ) + 3 5 *B ( 3 ) ** 2 + 10 5 *B ( 1 ) **4 ) *F 

ELSE DO 

IF (K . EQ. 9 )THEN DO 

C= ( - A**9 + 3 6 *A**7*B ( 1 ) - 8 4*A**6 *B ( 2 ) + 12 6 *A**5 *B ( 3 ) - 12 6 * 
*A**4*B ( 4 ) + 84*A**3*B ( 5 ) - 3 6*A**2*B ( 6 ) + 9*A*B ( 7 ) -B ( 8 ) - 
* 3 7 8 *a**5*B ( 1 ) **2 + 1230 *A**4 *B ( 1 ) *B ( 2 ) - 12 6 0 *A** 3 *B ( 1 ) * 
*B ( 3 ) + 7 5 6 *A**2*B ( 1 ) *B ( 4 ) - 2 5 2*A*B ( 1 ) *B ( 5 ) + 3 6 *B ( 1 ) *B ( 6‘) ) 

C= C+ 12 6 0*A**3 *B ( 1 ) **3 - 3 7 8 0*A**2*B ( 1 ) **2* 
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*B (2 ) - 840*A**3*B (2 ) '"*2+ 1890 *A* 

*B ( 1 ) **2*B ( 3 ) + 12 6 0*A**2*B ( 2 ) *B ( 3 ) - 3 7 8 * B ( 4 ) * 

*B ( 1 ) "'"2 - 404''' A "B ( 2 ) *B ( 4 ) + 84*B ( 2 ) *B ( 5 ) - 945 *A* 

" B ( 1 ) **4 + 2520 "A"B ( 1 ) "B ( 2 ) ** 2+1260 *B ( 1 ) ** 3 "B ( 2 ) - 3 1 5 * A* 
*B ( 3 ) **2 - 1260*B ( 1)*B (2 )*B (3 ) + 126*B ( 3 )*B (4) -280*B (2 ) ** 
C=C*F 



END IF 
END IF 
END IF 
END IF 
END IF 
END IF 
END IF 
END IF 
END IF 
END IF 
RETURN 
END 

$ ENTRY 



3 
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section Examples subsection All failure rates are dissimilar 



Shorthand notation : 

EXP ( . 3 ) + EXP ( .32)+ EXP ( . 4 ) + EXP ( . 6 ) + EXP (.62) 
Output 



MULTIPLICITY 

1 

1 

1 

1 

1 



TIME = 2.00000 

FAILURE RATE 
0.30000 
0.32000 
0.40000 
0.60000 
0.62000 

POLYNOMIAL COEFFICIENT FOR LAMDA 1 IS 
alO = 248.000000000000 

POLYNOMIAL COEFFICIENT FOR LAMDA 2 IS 
a20 = -332.142857142857 

POLYNOMIAL COEFFICIENT FOR LAMDA 3 IS 
a30 = 101.454545454545 

POLYNOMIAL COEFFICIENT FOR LAMDA 4 IS 
a40 = -70.857142857143 

POLYNOMIAL COEFFICIENT FOR LAMDA 5 IS 
a50 = 54.545454545454 



P (T > 2.00) 



0.998171698511 
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subsection All failure rates are identical : 



Shorthand notation : 

EXP (2.1) + ... + EXP (2.1) (10 times) 



Output 



TIME = 3.20000 

FAILURE RATE MULTIPLICITY 

2.10000 10 

POLYNOMIAL COEFFICIENT FOR LAMDA 1 ARE : 



alO 


= 


1.000000000000 


all 


= 


2.100000000000 


al2 


= 


2.205000000000 


al3 


= 


1.543500000000 


al4 


= 


0.810337500000 


al5 


- 


0.340341750000 


al6 


= 


0.119119612500 


al7 


= 


0.035735883750 


al8 


= 


0.009380669484 


al9 




0.002188822880 



P (T > 3.20) = 0.857717853017 



79 



subsection Some similar and some dissimilar failure rates 



Shorthand notation : 



EXP (6.5)+ 


... + 


EXP (6 . 5) 


(5 


times ) 


EXP (4.5)+ 


... + 


EXP (4. 5) 


(5 


times ) 


EXP (3.5)+ 


... + 


EXP (3.5) 


(5 


times ) 


EXP (2.5)+ 


... + 


EXP (2. 5) 


(5 


times ) 


EXP (1.5)+ 


... + 


EXP (1.5) 


(5 


times ) 


EXP ( .50)+ 


... + 


EXP ( .50) 


(5 


times ) 



Output 



TIME = 8.00000 



FAILURE RATE 

6.50000 

4.50000 

3 . 50000 

2.50000 

1.50000 
0.50000 



MULTIPLICITY 

5 

5 

5 

5 

5 

5 



POLYNOMIALS COEFFICIENTS FOR LAMDA 1 ARE : 
alO = -0.033655633163 

all = -0.015984638484 

al2 = -0.002968083834 

al3 = -0.000255685577 
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a!4 



-0.000008633539 



POLYNOMIAL 
a20 = 
a21 = 
a22 = 
a23 = 
a24 = 



COEFFICIENTS FOR LAMDA 2 ARE 
78824.407500664440 
28606.965700048750 
4242.082627828356 
308.298280647481 
9.469913057431 



POLYNOMIAL 
a30 = 
a31 = 
a32 = 
a33 = 
a34 = 



COEFFICIENTS FOR LAMDA 3 ARE 
-1939257.342698666000 
-833405.255639648600 
-199021.555234616700 
-18294.927183023870 
-1641.852439502142 



POLYNOMIAL 
a40 = 
a41 = 
a42 = 
a43 = 
a44 = 



COEFFICIENTS FOR LAMDA 4 ARE 
2591869.114320258000 
-558729.118591437400 
340359.641443397200 
-14083.275902489360 
4142.139971320400 



POLYNOMIAL 
a50 = 
a5 1 = 
a52 = 
a53 = 



COEFFICIENTS FOR LAMDA 5 ARE 
-736153.955891001200 
310180.360064797400 
-58447.530039485250 
5362.158719218832 
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a54 



-297.897706623269 



POLYNOMIAL 


COEFFICIENTS FOR LAMDA 6 ARE 


a60 = 


4718.810424379766 


a61 = 


-1542.762826749752 


a62 = 


203.749271368373 


a63 = 


-12.977288263887 


a64 = 


0.350737520646 



P(T > 8.00) = 0.998940885595 
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